Discrete charge patterns, Coulomb correlations and interactions in protein solutions 
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The effective Coulomb interaction between globular proteins is calculated as a function of mono- 
valent salt concentration c s , by explicit Molecular Dynamics simulations of pairs of model proteins 
in the presence of microscopic co and counterions. For discrete charge patterns of monovalent sites 
on the surface, the resulting osmotic virial coefficient B2 is found to be a strikingly non-monotonic 
function of c s . The non-monotonicity follows from a subtle Coulomb correlation effect which is 
completely missed by conventional non-linear Poisson-Boltzmann theory and explains various ex- 
perimental findings. 
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A more fundamental understanding of the interac- 
tions between nano-sized biomolecules is critical to the 
long-term advance of modern biomedical research [Q. 
The best strategy for a predictive calculation is to 
study simple coarse-grained models where effects can be 
clearly separated and approximations can be systemat- 
ically tested. While for micron-sized colloidal particles 
such coarse-grained models have led to a quantitative un- 
derstanding of the effective interactions g] , the challeng- 
ing question is how far this concept can be transferred to 
nano-particles. 

A particular issue is the aggregation and crystallization 
of globular proteins in solution, driven by their mutual 
interactions, including steric repulsion, van der Waals at- 
traction, Coulombic interactions, hydration forces, hy- 
drophobic attraction and depletion forces Most of 
these are effective interactions which depend sensitively 
on solution conditions. In particular Coulombic forces 
are functions of pH (which determines the total charge 
of the proteins) and of electrolyte concentration, which 
controls the Debye screening length A d , and hence the ef- 
fective range of Coulombic interactions. This dependence 
on solution conditions is exploited in "salting-out" exper- 
iments where large salt concentrations are used to trigger 
protein crystallization, a crucial step towards the deter- 
mination of their structure by X-ray diffraction || . While 
the forces acting between micro-sized colloidal particles 
are dominated by generic interactions, and are directly 
measurable by optical means @-^J, the interactions be- 
tween globular proteins are highly specific at short range, 
and are less directly accessible. One possible indirect de- 
termination of the total force between two proteins may 
be achieved via measurements of the osmotic equation of 
state by static light scattering, which in the low protein 
concentration regime yields the value of the second os- 
motic virial coefficient B2 J7j,|| . The variation of B2 with 
solution conditions yields valuable information on the un- 
derlying effective protein pair interactions. Moreover, it 
has been shown empirically that there is a strong corre- 
lation between the measured values of B2 and the range 
of solution conditions under which protein crystallization 
is achieved g-juj. 

This report focuses on the effective interactions be- 



tween globular proteins mediated by the microscopic co 
and counterions, and on the resulting T$i- The con- 
ventional Derjaguin-Landau-Verwey-Overbeek approach 
PH , borrowed from colloid science, leads one to expect 
that B2 will monotonically decrease as the concentration 
of salt increases, since higher salt concentrations lead to 
enhanced screening (i.e. reduction of Ad), and hence to 
a decrease of the effective protein diameter. This be- 
havior rests on the standard "coarse-grained" model of 
uniformly charged colloids and smoothed local densities 
of the microions. We show that the discrete nature of 
the protein surface charge distribution, together with the 
Coulomb correlations between all charges involved, lead 
to a striking non-monotonic variation of B2 with salt con- 
centration c s . The occurrence of a minimum of B2 as a 
function of c s has recently been reported in lysozyme so- 
lutions for c s — 0.3 M |T^] and in Apoferritin solutions for 
c s = 0.15 M (ll|. Related experimental findings are non- 
monotonic variations of other quantities which strongly 
correlate with B2 0,|lj] such as the interaction param- 



,the cloud point temperature |l7],[l8| , and the 
1S[|. All these trends can be qualitatively un- 



eter [[15|,|16 
solubility 
derstood by our calculation. 




FIG. 1. Snapshot of a typical MD-generated microion con- 
figuration around two proteins, separated by r — 1.7a p . The 
proteins carry 15 discrete charges e; monovalent salt molar- 
ity is c s = 0.206Mol/l. The globular protein molecules are 
shown as two large gray spheres. The embedded small dark 
spheres on their surface mimic the discrete protein charges 
in the DCM model. The small gray spheres are counterions, 
while the black spheres are coions. 
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We consider two spherical proteins of diameter a p , each 
carrying a total charge Ze (where Z depends on pH), sur- 
rounded by monovalent co and counterions, assumed to 
have identical diameters a c . The solvent (water) is as- 
sumed to be a dielectric continuum of permittivity e; this 
simplification, which ignores the molecular granularity of 
the solvent, amounts to the standard "primitive" model 
of ionic solutions [^o) . 

In the case of highly charged colloidal particles, the 
total charge Ze is usually assumed to be uniformly dis- 
tributed on the surface, a situation which will be referred 
to as the "smeared charge model" (SCM). This simpli- 
fication is much less justified for the smaller, weakly 
charged proteins (where Z ~ 10). We have hence 
adopted a second, discrete charge model (DCM) where 
Z monovalent discrete point charges are distributed over 
the surface of a sphere of diameter a p — 2A = 0.96er p 
(i.e. slightly inside the surface of the protein), in such a 
way as to minimize the electrostatic energy of the dis- 
tribution; the resulting pattern does not correspond to 
the real charge distribution on any specific protein, but 
does provide a well-defined discrete model for comparison 
with the SCM, and between different values of Z. At this 
stage the two models (SCM and DCM) involve only ex- 
cluded volume and bare Coulombic interactions (reduced 
by a factor 1/e to account for the solvent) between all 
particles (proteins and microions). However, in view of 
the large size asymmetry, the effective force between the 
proteins, which ultimately determines the second virial 
coefficient, involves a statistical average over microion 
configurations in the field of two fixed proteins j2^] . For 
distances r > a p between the centres of proteins 1 and 2, 
the total force F\ = — F 2 acting on each of the proteins is 
the sum of three contributions, F\ = F± + F± + F± , 
where is the direct Coulomb repulsion between the 

~~*(2) 

charges on the two proteins, F± is the microion induced 

electrostatic force and F± is the depletion force due to 
the imbalance of the microion osmotic pressure acting on 
opposite sides of the proteins [ pl[ . Both F^ and F^ 
are averages over microion configurations; according to 
the contact theorem F± is directly related to the inte- 
gral of the microion contact density over the surface of 
the protein p2jp3]. The statistical averages leading to 



Fi and F^' were computed using Molecular Dynamics 
(MD) simulations. The two proteins were placed sym- 
metrically with respect to the centre along the body di- 
agonal of a cubic simulation cell of length L = 4a p , which 
also contained monovalent co and counterions in numbers 
determined by their bulk concentrations; periodic bound- 
ary conditions were adopted. The choice of L was made 
to ensure that the box length is much larger than the 
range of the total (effective) protein-protein interaction, 
so that the results would be independent of L. For our 
model to be a rough representation of lysozyme, we chose 
(T p = Anm, and Z = 6, 10 and 15, corresponding to three 
different values of the solution pH. The microion diame- 



r(3 



ter is <j c = 0.267nm. Note that the SCM always implies 
vanishing multipole moments, whereas within the present 
DCM, the only charge pattern with a non-vanishing 
dipole moment is that for Z = 15. A snapshot of a typical 
equilibrium microion configuration around two proteins 
is shown in Fig.Q, for the case Z = 15. Note that the di- 
mensionless Coulomb coupling parameter for a protein- 
counterion contact, namely T = e 2 /[e/cBT(A + cr c /2)] for 
the DCM, and T = 2Ze 2 /[ek B T(a p + a c )] for the SCM, 
are comparable and of the order of Y — 3 at room tem- 
perature. The total force F\ = —F2 depends only on 
the centre-to-centre distance r for the SCM, but is also 
a function of the orientations of the two charge patterns 
of the DCM, embodied by two unit vectors uj\ and Q 2 \ 
F\ = F\(r, Wi, (2%). The anisotropy of the force turns out 
to be relatively weak. The effective radial pair interac- 
tions between proteins, V(r), follow from integration of 
the radial projection of an orientationally averaged force 
F\ along the centre-to-centre vector r, according to: 



V(r) 



dr'(—7 ■ F(r',uji,uj 2 )} Ul u 



(1) 



Here (. . ■}u, 1 ,u 2 refers to a statistical average over mutual 
orientations of the two proteins J24J. The second virial 
coefficient in units of its value 2tto^/3 for hard spheres 

of diameter a p , B 2 = B 2 / B 2 , can then be proven to 
be given by: 



B* 
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exp{-V(r)/k B T}}, (2) 



a result formally identical to that valid for spherically 
symmetric forces. Results for B 2 as a function of salt 
concentration are shown in Fig.^ for the SCM and DCM 
models, with three values of the total protein charge. In 
order to obtain values of B2 comparable to measured 
virial coefficients, we have taken short-range attractions 
between proteins into account, by adding to the effective 
Coulomb potential in Eq. (g) a "sticky" hard-sphere po- 
tential of the Baxter form pa ] , with potential parameters 
5 = 0.02cr p and r = 0.12, which arc known to yield rea- 
sonable osmotic data for lysozyme solutions ||^6| in the 
high salt concentration regime, where Coulombic inter- 
actions are essentially screened out. 

The key result, illustrated in Fig.||, lies in the consider- 
able qualitative difference between the predictions of the 
SCM and the DCM models for the variation of B% with 
monovalent salt concentration c s , irrespective of the total 
protein charge Ze. While the SCM (dashed curves) pre- 
dicts a monotonic decay of B 2 with c s , the DCM leads 
to a markedly non-monotonic variation, involving an ini- 
tial decay towards a minimum followed by a subsequent 
increase to a maximum and a final decrease towards a 
high c s value similar to that predicted by the SCM. The 
location of the minimum and of the maximum shift to 
higher values of c s for larger protein charges Z . 

The origin of the non-monotonic variation of B 2 with 
c s can be traced back to the dependence of the effective 
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(screened) Coulomb interaction on salt concentration as 
shown in the inset of Fie .| for Z = 10. While the spher- 
ically averaged, repulsive effective potential V(r) of the 
DCM is initially strongly reduced as c s is increased, its 
amplitude and range increase very significantly at inter- 
mediate concentrations (c s ~ 1M/1), before it nearly van- 
ishes at the highest salt concentrations. Note that V(r) 
becomes even slightly attractive at contact (r = <j p ) for 
c s ~ 2M/1. The enhanced effective Coulomb repulsion 
at intermediate salt concentrations cannot be rational- 
ized in terms of simple mean-field screening arguments; 
it is caused by a subtle correlation effect which leads to 
the non-monotonic behavior of B2 within the DCM. The 
protein-microion correlations are of a sufficiently different 
nature in the SCM, to lead to a much more conventional, 
monotonic decay of B2 with c s , similar to that expected 
from a linear screening picture. 




FIG. 2. Nor- 

malized second virial coefficient = B2/B^ s of a protein 
solution versus added salt mofarity. Results are shown for 
protein charges Z = 6 (dashed lines) , Z = 10 (solid fines) and 
Z = 15 (dot-dashed fines). The fines with (without) symbois 
correspond to the SCM (DCM) modef. The inset shows the 
effective protein-protein interaction V(r) in the DCM model 
versus separation distance r for Z — 10. Various symbois in 
the inset reiate to the different added saft mofarities, indicated 
in the iegend. 

Even though the effective Coulomb potential between 
proteins is of small amplitude, only a few percent of the 
thermal energy fc^T, the effect on B2 is dramatically en- 
hanced by the presence of the strong short-range attrac- 
tive component due to van der Waals and hydrophobic 
interactions, which we have included in the form of the 
Baxter "sticky" sphere potential. This potential is in- 
dependent of salt-concentration, and has no influence on 
the qualitative dependence of B2 on c s . 

In order to gain further insight into the physical mecha- 
nism responsible for the unusual variation of the effective 
Coulomb potential and of B2 with salt concentration, we 
have investigated in detail the local microion density in 
the immediate vicinity of the protein surface. The radial 



microion density profile p(r) — p+(r) + p~(r) around a 
single isolated protein is shown in FigJ|, for Z = 10, and 
two salt concentrations (the profiles are orientationally 
averaged in the case of the DCM). At the lower salt con- 
centration (c s = 0.206M/0 the SCM and DCM models 
both yield an accumulation of the microion density near 
contact, in semi-quantitative agreement with the predic- 
tion of standard Poisson-Boltzmann (PB) theory. At the 
higher salt concentration, however, there is a marked de- 
pletion in the microion density, signaled by a minimum 
of p(r) well below the asymptotic bulk value. This cor- 
relation effect is of course absent in the (non-linear) PB 
theory, which always predicts a monotonically decreasing 
density profile p(r). Note however a significant difference 
between the SCM and DCM profiles. While the latter 
predicts a contact value p c (r — (a p + a c )/2) larger than 
the bulk value, SCM predicts a much stronger microion 
depletion near contact. This finding illustrates the sensi- 
tivity of correlation effects to the assumed charge pattern 
at the surface of a protein: taking into account the dis- 
creteness of the surface charges leads to a significant re- 
duction of microion depletion at contact, compared to the 
simplified picture of a uniformly smeared charge (SCM) . 
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FIG. 3. Total density profiies p(r) = p+(r) — p-(r) of smalf 
ions around a singie protein, for salt molarities c s — 0.206 
(bottom set of curves) and c s = 0.824 (upper set of curves). 
The solid and dashed lines are simulation results for DCM 
and SCM models respectively, while the dot-dashed lines are 
predictions of non-linear Poisson-Boltzmann theory. 

Next, consider the influence of a second near-by protein 
on the microion distribution near contact. We have com- 
puted the difference between "inner" and "outer" shell 
microion contact densities, as schematically illustrated in 
the inset to Fig.||. The local microion density is no longer 
spherically symmetric, due to the interference of the elec- 
tric double-layers associated with the two proteins. The 
difference Ap = pi„ — p out between the mean number of 
microions within a fraction of a spherical shell of radius 
R = 0.6(Tp subtended by opposite 60° cones is plotted in 
Fig.^ versus salt concentration. Ap is always positive, in- 
dicating that microions (in fact mostly counterions) tend 
to cluster in the region between the proteins, rather than 
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on the opposite sides, as one might expect due to the 
enhanced lowering of the electrostatic energy for counte- 
rions shared between the two proteins. However, there is 
a very significant difference in the variation of Ap with c s , 
between the SCM and the DCM models. Both exhibit 
similar behavior for c s < 0.5M/1, with a small maxi- 
mum around 0.2M /I. Beyond 0.5M /I, however, the SCM 
predicts a monotonic decrease of Ap, while the DCM 
leads to a sharp peak in Ap for c s ~ XM/l. This highly 
non-monotonic behavior clearly correlates with the non- 
monotonicity of B 2 evident from Fig.^. The excess num- 
ber of microions between the two proteins leads to an 
imbalance in osmotic pressure, which is the origin of the 
increased repulsion between proteins around c s = 1M/1, 
as shown in the inset of FigJ^. For stronger Coulomb cou- 
pling, as is the case for highly charged colloidal particles 
in the absence of salt, the above depletion mechanism 
is inverted, and leads to a depletion attraction between 
the particles plfl , rather than to the enhanced repulsion 
found here in the case of relatively weakly charged pro- 
teins. 
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FIG. 4. Microion density imbalance Ap versus salt mo- 
larity for protein charge Z = 10 and separation r = 1.2a p . 
The solid and dashed lines correspond to the DCM and SCM 
models respectively. The inset shows the angular range over 
which Ap is averaged (see text). 

The main finding of the present work is that the sec- 
ond osmotic virial coefficient of protein solutions has a 
non-monotonic dependence on salt concentration if the 
charge pattern on the protein surface is discrete (as is 
the case for real proteins) rather than uniformly smeared 
out, as usually assumed in the related case of charge- 
stabilized colloidal dispersions, involving much larger 
particles. The lesson to be learned from this finding 
is that one must be cautious in attempting to extend 
coarse graining concepts and approximations, developed 
and routinely used on the colloidal scale, to the nano- 
metric scale of proteins. The discreteness of the charge 
pattern is crucial to obtain non-monotonic behavior of 
B2, which is a subtle Coulomb correlation effect, totally 



missed by non-linear PB theory. 

We chose our simple models to help highlight and sep- 
arate the effects of discrete charge patterns and Coulomb 
correlations. Extending our MD calculations to the more 
complex (pH dependent) charge patterns of realistic pro- 
teins |27j is technically straightforward. We expect that 
the physical mechanism leading to enhanced protein re- 
pulsion at intermediate salt concentration, which is illus- 
trated by the microion density imbalance shown in FigJ^, 
will carry over. Since the second osmotic virial coefficient 
determines much of the excess (non-ideal) part of the 
chemical potential of semi-dilute protein solutions, it is 
anticipated that the non-monotonicity of B 2 may have a 
significant influence on protein crystallization from such 
solutions in the course of a "salting-out" process. The 
non-monotonic behavior also suggests the possibility of 
an inverse, "salting-in" effect, whereby a reduction of salt 
concentration may bring B2 into the "crystallization slot" 
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